FFT,NTT

定义

多项式

  1. 系数表示:
      设 $A(x)$ 表示一个 $n-1$ 次多项式,所有项的系数组成的n维向量 $( a_{0},a_{1},a_{2},…,a_{n-1})$ 唯一确定了这个多项式,即
  2. 点值表示:
      将一组互不相同的 $n$ 个 $x$ 带入多项式,得到 $n$ 个值设它们组成的 $n$ 维向量分别为 $( x_{0}, x_{1}, x_{2},…, x_{n-1} )$ $( y_{0}, y_{1}, y_{2},…, y_{n-1} )$ ,即
  3. 求值与插值:
      首先我们知道一个 $n-1$ 次多项式在 $n$ 个不同点的取值唯一确定了多项式。朴素按照定义求点值复杂度 $O(n^{2})$,已知点值表示求系数表示,朴素插值法也为 $O(n^{2})$。
  4. 多项式乘法
      我们定义两个多项式 相乘的结果为 $C(x)$,则 直接相乘复杂度 $O(n^{2})$,如果一个个去算的话,要花费 $O(n^{2})$ 的时间才可以完成,但是,这是在系数表示下计算的,如果转换成点值表示,知道了 $A(x), B(x)$ 的点值表示后,由于只有 $n$ 个点,就可以直接将其相乘,在 $O(n)$ 的时间内得到 $C(x)$ 的点值表示。

单位根(n为2的正整数次幂)

  在复平面上,以原点为圆心,$1$ 为半径作圆(单位圆),以原点为起点,单位圆的 $n$ 等分点为终点,作 $n$ 个向量,设所得幅角为正且最小的向量对应的复数为 $\omega_{n}$ ( $n$ 次单位根),由乘法的定义(模长相乘,幅角相加)可知,$n-1$ 个向量分别对应 $\omega_{n}^{2}, \omega_{n}^{3} ,…,\omega_{n}^{n}$,其中$\omega_{n}^{n}=\omega_{n}^{0}=1$。
性质:

快速傅里叶变换

  考虑多项式 $A(x)$ 的表示,将 $n$ 次单位根的 $0$ 到 $n-1$ 次幂带入多项式的系数表示,所得点值向量 $(A(\omega_{n}^{0}),A(\omega_{n}^{1}),…,A(\omega_{n}^{n-1}))$ 称为其系数向量 $(a_{0},a_{1},a_{2},…,a_{n-1})$ 的离散傅里叶变换。朴素求法复杂度仍为$O(n^{2})$,所以寻求更高效的方法。
考虑将多项式按照系数下标的奇偶分为两部分

假设 $k<\frac{n}{2}$,现在求 $A(\omega_{n}^{k})$

对于

即当 $k<\frac{n}{2}$ 时,

注意到,当 $k$ 取遍 $[0,\frac{n}{2}-1]$ 时, $k$ 与 $k+\frac{n}{2}$ 取遍 $[0,n-1]$。
也就是说,如果已知 $A_{1}(x),A_{2}(x)$ 在 $\omega_{\frac{n}{2}}^{0},\omega_{\frac{n}{2}}^{1},…,\omega_{\frac{n}{2}}^{\frac{n}{2}-1}$ 处的点值,就可以在 $O(n)$ 的时间内得到 $A(x)$ 在 $\omega_{n}^{0},\omega_{n}^{1},…,\omega_{n}^{n-1}$ 处的点值。
而 $A_{1}(x),A_{2}(x)$ 的问题相当于原问题的一半,分治边界为一个常数项,复杂度为 $O(nlogn)$ ,这就是最常用的算法 $Cooley-Tukey$ 算法,即蝶形算法。

傅里叶逆变换

  将点值表示的多项式转换为系数表示,同样可以使用快速傅里叶变换,这过程叫傅里叶逆变换。
  设 $( y_{0},y_{1},y_{2},…,y_{n-1} )$ 为 $( a_{0},a_{1},a_{2},…,a_{n-1} )$ 的傅里叶变换。
  考虑另一个向量 $( c_{0},c_{1},c_{2},…,c_{n-1} )$ ,满足

即多项式

在 $\omega_{n}^{0},\omega_{n}^{-1},…,\omega_{n}^{-(n-1)}$ 处的点值表示。
将上式展开

考虑一个式子

当 $k\neq0$ 时,两边同乘 $\omega_{n}^{k}$ ,即

两式相减得

当 $k=0$ 时,$S(\omega_{n}^{k})=n$
继续考虑

当 $j=k$ 时, $S(\omega_{n}^{j-k})=n$ ,否则 $S(\omega_{n}^{j-k})=0$ ,即$c_{i}=na_{i}$ ,所以

综上,将单位根的倒数代替单位根,做一次FFT,再将结果除以 $n$ ,就是逆变换的结果。

实现(迭代法)

  1. 二进制翻转:考虑FFT分治到边界时,每个数的顺序及其二进制位。
    fft
    即分治后的下标等于原下标的二进制位翻转,代码实现:

    1
    2
    3
    4
    void init(int n){
    int k=0; while((1<<k)<n)k++;
    for(int i=1;i<n;i++)rev[i]=rev[i>>1]>>1|((i&1)<<(k-1));
    }
  2. 蝴蝶操作:考虑合并两个子问题的过程,假设 $A_{1}(\omega_{\frac{n}{2}}^{k})$ 和 $A_{2}(\omega_{\frac{n}{2}}^{k})$ 保存在 $a(k),a(\frac{n}{2}+k)$ 中,令 $t=\omega_{n}^{k}*a(\frac{n}{2}+k)$ ,则

    代码

    1
    2
    3
    4
    5
    6
    7
    8
    9
    10
    11
    12
    13
    14
    15
    16
    17
    18
    19
    20
    21
    22
    23
    24
    25
    26
    27
    28
    29
    30
    31
    32
    33
    34
    struct cp{
    double x,y;
    cp(double x=0,double y=0):x(x),y(y){}
    cp operator +(const cp &d)const{return cp(x+d.x,y+d.y);}
    cp operator -(const cp &d)const{return cp(x-d.x,y-d.y);}
    cp operator *(const cp &d)const{return cp(x*d.x-y*d.y,x*d.y+y*d.x);}
    cp operator /(const double &d)const{return cp(x/d,y/d);}
    };
    struct FFT{
    int rev[maxn];
    void init(int n){
    int k=0;while((1<<k)<n)k++;
    for(int i=1;i<n;i++)rev[i]=rev[i>>1]>>1|((i&1)<<(k-1));
    }
    void trans(cp *a,int n,int f){
    for(int i=0;i<n;i++)if(i<rev[i])swap(a[i],a[rev[i]]);
    for(int l=2;l<=n;l<<=1){
    int m=l/2; cp omega=cp(cos(2*pi/l),sin(2*pi/l)*f);
    for(cp *p=a;p!=a+n;p+=l){
    cp w=cp(1,0);
    for(int i=0;i<m;i++){
    cp t=p[i+m]*w;
    p[i+m]=p[i]-t; p[i]=p[i]+t;
    w=w*omega;
    }
    }
    }
    }
    void dft(cp *a,int n){trans(a,n,1);}
    void idft(cp *a,int n){
    trans(a,n,-1);
    for(int i=0;i<n;i++)a[i]=a[i]/n;
    }
    }fft;
0%